Efficient construction of maximally localized photonic Wannier functions: 
locality criterion and initial conditions 
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Wannier function expansions are well suited for the description of photonic-crystal-based defect 
structures, but constructing maximally localized Wannier functions by optimizing the phase degree 
of freedom of the Bloch modes is crucial for the efficiency of the approach. We systematically analyze 
different locality criteria for maximally localized Wannier functions in two-dimensional square and 
triangular lattice photonic crystals, employing (local) conjugate- gradient as well as (global) genetic- 
algorithm-based, stochastic methods. Besides the commonly used second moment (SM) locality 
measure, we introduce a new locality measure, namely the integrated modulus (IM) of the Wannier 
function. We show numerically that, in contrast to the SM criterion, the IM criterion leads to 
an optimization problem with a single extremum, thus allowing for fast and efficient construction 
of maximally localized Wannier functions using local optimization techniques. We also present an 
analytical formula for the initial choice of Bloch phases, which under certain conditions represents 
the global maximum of the IM criterion and, thus, further increases the optimization efficiency in 
the general case. 
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I. INTRODUCTION 

Photonic crystals (PhCs) remain to attract a consider- 
able attention of the scientific community due to their 
unique properties and potential technological applica- 
tions To a large extent, PhCs applications are based 
on the photonic bandgap effect and involve sophisticated 
defect structures, cavities and wave guides in the periodic 
crystal host. The full dynamics of the electromagnetic 
field in such structures may be studied in principle rigor- 
ously by direct numerical solution of Maxwell's equations 
using the finite difference time domain (FDTD) method 
For the calculation of stationary modes, however, the 
numerical effort may be substantially reduced by using 
the Galerkin method [3], i.e., by expanding the electro- 
magnetic field in terms of an appropriate orthonormal set 
of basis functions which renders the stationary electro- 
magnetic wave equation as a discrete matrix eigenvalue 
problem. In this case the proper choice of the set of basis 
functions is crucial in order to obtain an accurate de- 
scription while keeping the dimension of the eigenvalue 
problem minimal. 

While for extended wave problems an expansion in 
terms of Bloch functions, the eigenmodes of the unper- 
turbed PhC, is natural [j-l6j, for the description of de- 
fect structures the use of Wannier functions as a ba- 
sis set is in principle superior, because Wannier 
functions may be constructed as being localized in space 
and are still an exact representation of the point symme- 
try group of the host PhC. Effective solutions developed 
for the electronic case, namely, maximally localized gen- 
eralized Wannier functions [ll|, [l2|, have been applied 
only recently to the electromagnetic case [HI, [H| . Since 
then, the theory of photonic Wannier functions has been 
applied to the analysis of 2D PhC cavities, waveguides 
[3, [H| , waveguide crossings [l6| and PhC heterostruc- 



tures [17]. The generalization of the approach to the 
case of 2D slab PhCs and 3D PhCs has been also re- 
ported in [l8|, EH- However, the construction of well- 
localized Wannier functions involves a multidimensional 
optimization problem [20|, with the arbitrary phase of 
each Bloch mode (defined below) as optimizations pa- 
rameters. Therefore, the practical importance of this ap- 
proach has been fairly limited up to now. In the present 
work we systematically analyze this optimization prob- 
lem employing local as well as global optimization pro- 
cedures. We propose a novel locality measure for the 
Wannier functions which allows for a highly efficient op- 
timization of the locality of the Wannier function. 

Wannier functions are defined as the Fourier transform 
of the Bloch modes 



B nk (r) = e , *» k e lkr u nk (r) 
with respect to the wave vector k, 



W„ R (r) = ^=J2 e" ikR B„ k (r). 



(i) 



(2) 



The Bloch mode is a solution of the corresponding wave 
equation in a periodic medium, where u n k(r) is a lattice- 
periodic envelope function, and we have explicitly de- 
noted the arbitrary phase n k of the Bloch function, the 
Bloch phase, n is the band index. For simplicity, in 
the present work we do not consider mixing of different 
bands [l6|, [2l| in the construction of Wannier functions 
from Bloch modes. 

From the definitions (pQ) , (j2j) it is seen that the Wannier 
functions would be (5(r)-localized in space only if the en- 
velopes u n k(r) were constant. By contrast, for arbitrary 
Bloch phases </> n k the Wannier functions have actually a 
rather large spatial extension due to the oscillatory char- 
acter of the u n k(r), a feature which becomes more and 
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more pronounced for higher band indices n. However, 
the gauge freedom in the Bloch phases may be employed 
to adjust, for each k, the value of n k so as to construct 
maximally localized Wannier functions with respect to 
some locality criterion. This constitutes a complex, mul- 
tidimensional optimization problem. The choice of the 
locality measure is not unique. The common choice used 
in the literature [13|, [2l|, |22| is to minimize the second mo- 
ment (SM) of the modulus square of the Wannier func- 
tion [ill fl2j|. Unfortunately, it turned out that the SM, 
as a functional of the set of Bloch phases n k, has multi- 
ple local minima, so that local optimization procedures, 
like the conjugate gradient method, tend not to find the 
global minimum and the locality of the Wannier function 
optimized in this way depends sensitively on the initial 
set of Bloch phases. Global optimization procedures are, 
on the other hand, slow and computationally exceedingly 
expensive. 

The purpose of this paper is a systematic analysis of 
different locality measures using both local and global op- 
timization methods. In Section [II] we construct SM opti- 
mized maximally localized Wannier functions, comparing 
conjugate gradient (local) and genetic- algorithm-based 
(global) methods. This analysis leads us in Section Hill to 
propose a new locality measure, namely the integrated 
modulus square (IM) measure, resulting in an optimiza- 
tion problem with a single extremum only. This allows 
one to use fast and efficient local optimization techniques 
to construct maximally localized Wannier functions. In 
Section [TV] we show that, if the Bloch modes conform 
certain conditions, it is possible to find an optimal set of 
Bloch phases with respect to the IM locality measure an- 
alytically. Although the required constrains will not be 
fulfilled in the general case, the Wannier functions calcu- 
lated using such an analytical set of phases show strong 
tendency towards localizations and can be used as an ef- 
ficient starting point for the numerical optimization. In 
Section [V] using several examples of PhC defect struc- 
tures, the quality of the constructed Wannier functions 
is demonstrated. Section [VTl concludes the paper. 



II. SECOND MOMENT OPTIMIZATION 

A. Definitions 

In what follows, we limit ourselves to the two- 
dimensional (2D) case, which is characterized by the pe- 
riodic dielectric functions e(r) = e(r + R), VR G £, with 
r = (x, y) denoting a 2D vector in the x-y-plane and 
R being a lattice vector of some 2D lattice C. In this 
case, the wave equation for time harmonic TM (trans- 
verse magnetic), E(r,t) = e~ luJt E(r), and TE (transverse 
electric), H(r,t) = e~ lujt H(r), polarization reads, respec- 



tively, 

£ E E(r) = --L\7 2 E(r) = ^E(v) (3) 
C H H(r) = -V-^ViT(r) = ^tf(r). (4) 

The wave operators Ce and Ch are hermitian with re- 
spect to the corresponding inner products: 

(f\9)E = I d\t{r)e{v)g{v) (5) 
Jv 

(f\9) H = f d 2 rr(r)g(r). (6) 
Jv 

V is the 2D volume of the crystal. It is important to 
mention, that the completeness and orthogonality of the 
Bloch modes i? n k(r) and H n k(r) translates into the com- 
pleteness and orthogonality of Wannier functions W n n(r) 
with respect to the corresponding inner product: 

(W n Ti\W n >ii'}E/H = <W#RR'- (7) 

Moreover, the translation property of the Wannier func- 
tions, 

WnR(r) = W n0 (r-R), (8) 

follows from the periodicity of the envelope functions 
^nk(r). 

In this work four different two-dimensional photonic 
crystals are analyzed for both fundamental polarizations. 
Square (Sq) and triangular (Tr) lattices of dielectric rods 
in air (D) and air rods in dielectric (A) are considered. In 
what follows we will refer to these systems as Sq-D, Tr- 
D, Sq-A and Tr-A, respectively. The radius of rods and 
the dielectric constant of dielectric material and air are 
chosen to be ro/a = 0.2, e = 12 and e = 1, respectively. 
We adopt the following definition of the second moment 
of the Wannier function W n n(r): 

5 n ({0nk}) = <Wn R |(r - T ) 2 \W n n) E / H , (9) 

with i*o being the Wannier center. Two positions of the 
Wannier center will be analyzed further, (i) in the cen- 
ter of the scatterer ("on-site") and (ii) in the geometri- 
cal center between four (three) scatterers in the case of 
square (triangular) lattices ("between"). For example, 
in a square lattice with lattice constant a, the Wannier 
center can be set to ro = R (on-site) or ro = R + T 
(between), where T = (0.5a, 0.5a). 

B. Conjugate gradient method 

The commonly adopted way of finding maximally lo- 
calized Wannier functions is to minimize the correspond- 
ing second moment S n with respect to the set of Bloch 
phases </> n k in the nth band. Note that, by definition 
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FIG. 1: Second moments S n (inverse locality) of the Wannier 
functions in the nth band, minimized by using the conjugate 
gradient method. Four different, randomly chosen initial sets 
of Bloch phases, A, B, C, D, were used for the CG optimiza- 
tion. Top: Sq-D crystal, TM polarization. Bottom: Tr-D 
crystal, TE polarization. 



of the SM, the regions far away from the Wannier cen- 
ter (where the Wannier functions usually have an in- 
creasingly complex structure) contribute to the SM with 
quadratically increasing weight. For this optimization we 
first apply a standard conjugate gradient (CG) method. 
This method tends to get trapped in the local extrema 
and, hence, requires a careful choice of the initial set of 
Bloch phases [21]. In figure [1] representative examples of 
the locality of the SM-optimized Wannier functions are 
shown for Sq-D (TM polarization) and Tr-D (TE polar- 
ization) structures. Here and in the following, the Bloch 
modes have been calculated using the plane wave expan- 
sion method [23]. The minimized SM (i.e., the inverse 
locality) of the Wannier functions in the first eight bands 
are displayed for four different, random distributions of 
the initial Bloch phases. As expected, the optimal local- 
ity obtained using the CG method depends crucially on 
the choice of the initial phases. This is an indication that 
the SM of the Wannier functions possesses several local 
minima. 







- V 


i + * T 







IM=0.25 
1/SM=2.23 



r, C O 



20000 40000 
generations 



IM=0.24 
1/SM=1.82 



40000 80000 120000 
generations 



0.3 
0.2 

0.1 


0.3 
0.2 
0.1 





FIG. 2: Locality (1/S n ) of the SM-optimized Wannier func- 
tions as a function of GA generation for the Sq-D structure, 
TM polarization, 3rd band (top) and for the Sq-A structure, 
TE polarization, 5th band (bottom). The dashed, green line 
shows the locality of the best- localized Wannier function in 
each generation. Every 100 generations, these Wannier func- 
tions served as a starting point for the subsequent CG opti- 
mization step (red crosses). On the right hand side, the mod- 
ulus square of the SM-optimized Wannier function is shown 
for an early (top: 40000th, bottom: 38300th) and a later (top: 
60000th, bottom: 100000th) generation, respectively. 



C. Genetic algorithm 

The proper choice of the initial set of Bloch phases is 
in general not trivial [21]. Even if such a choice can be 
justified, one can never be sure that the resulting SM- 
optimized Wannier functions indeed correspond to the 
global minimum of the SM. To solve the global opti- 
mization problem and to examine the validity of the SM 
optimization in more detail, we have applied a global, 
stochastic-based optimization technique, namely a ge- 
netic algorithm (GA) method [24] . Taking the biologic 
evolution in nature as a model, the GA method works 
with a population of individuals which pass through a 
selection procedure and can reproduce themselves. Each 
Wannier function represents an individual. The set of 
Bloch phases, which determines the Wannier function, is 
represented as a large, binary string. The GA method 
starts with a population of random Wannier functions 
and passes them through a selection procedure where 



only that one half of the Wannier functions are retained 
which are most strongly localized with respect to the 
given locality criterion ("most fit individuals). In a sec- 
ond step, these survived individuals are allowed to re- 
produce themselves by randomly mixing their strings of 
phases, thus passing their attributes to the offsprings. 
Together with their off-springs, the survived individuals, 
which correspond to better localized Wannier functions, 
comprise the new generation. By iterating this proce- 
dure over several thousands of generations the algorithm 
will converge slowly but definitely towards the global ex- 
ternum. Once the GA procedure has reached the valley 
of the global extremum, the CG method should be ap- 
plied subsequently to the GA algorithm in order to ac- 
celerate the convergence and improve the accuracy of the 
solution [24j . 

In figure [2] the evolution of the G A results is depicted 
for two representative systems and polarizations. Every 
100 generations a Wannier function with highest locality 



4 




5 









*-« Phases A 










°- -° Phases B 










«™» Phases C 










Phases D 






I 


i i i i 


i 


1 


2 


3 4 5 6 


7 E 



band index 



FIG. 3: Locality (T n ) of the Wannier functions optimized with 
respect to the IM locality measure using the CG method. 
Four different randomly chosen initial sets of Bloch phases 
were used for the CG optimization. Top: Sq-D crystal, TM 
polarization. Bottom: Tr-D crystal, TE polarization. 



in the current population was taken as a starting point 
for the subsequent CG optimization. Over the first sev- 
eral thousand generations the locality of the resulting 
Wannier functions is varying strongly, indicating hopping 
of the solution among different local minima due to the 
stochastic nature of the algorithm. At the top panel of 
figure [2] one can observe how the algorithm is stuck in 
a local minimum over several thousands of generations, 
before it escapes and reaches the global minimum valley 
at around the 50000th generation. An improvement of 
the locality for later generations is clearly seen from the 
modulus square of the optimized Wannier functions (fig- 
ure [2j right panels). The discontinuous nature of the GA 
method ensures with stochastic certainty that the global 
minimum of the SM is found, providing the best local- 
ization of the Wannier functions with respect to a given 
locality measure. At the same time, however, the numer- 
ical load of the GA method exceeds the one of the CG 
method by far, making it inappropriate for routine appli- 
cation for an efficient construction of maximally localized 
Wannier functions. 



III. INTEGRATED MODULUS OPTIMIZATION 

The complicated structure of the Wannier functions at 
large distances, which is expressed by several local min- 
ima of their SMs and the associated difficulties in the 
construction of maximally localized Wannier functions, 
motivates the search for a simpler criterion for the lo- 
cality of Wannier functions. Here, we introduce a new 
criterion, the integrated modulus square (IM), defined as 

X n {{M) = / d 2 rW; R (r)I(r)W nR (r), (10) 

jug 

where the integration region is the first unit cell (UC) 
around the Wannier center (UC). We choose the function 
X(r) in such a way, that the IM is equal to unity for 
a Wannier function which is completely confined within 



such a unit cell, i.e., 

X(r) "\l forTE • 

A well localized Wannier function corresponds to a large 
IM, and one needs to maximize the IM in order to obtain 
the maximally localized Wannier functions. The IM is 
a very sharp criterion, since it does not depend on the 
structure of the Wannier functions outside of the integra- 
tion region at all. 

We examined the IM locality measure as a localization 
criterion for the same four different physical systems and 
both polarizations as it was done in the SM case. To that 
end we construct maximally localized Wannier functions 
using the same four different, randomly chosen initial sets 
of Bloch phases as in Section [TTJ Representative examples 
of CG optimization are shown in figure [3] In contrast 
to the SM case, the localities of the resulting Wannier 
functions coincide for all four sets of phases. This is the 
case for all considered systems and polarizations, strongly 
indicating that the IM locality measure does not possess 
any local extrema. 

To support this hypothesis, we applied the GA method 
to solve the optimization problem. For all considered 
structures no significant variation of the locality has been 
observed for different GA generations. Representative 
examples are shown in figure SJ One can clearly see, that 
over many thousands of generations the locality of the 
Wannier functions optimized with respect to IM criterion 
stays constant. This proves numerically that the consid- 
ered optimization problem possesses a single extremum, 
making the procedure independent on the choice of the 
initial set of Bloch phases. As a consequence, the use of 
the IM as a locality criterion together with the CG as 
an optimization method represents a fast as well as reli- 
able method for the construction of maximally localized 
Wannier functions. Figures [5] and [6] show the maximally 
localized Wannier functions for the first several bands of 
the Sq-D (TM polarization) and Sq-A (TE polarization) 
structures, respectively. In both cases the IM criterion 
was used along with the CG method. The optimized 
Wannier functions demonstrate good locality which de- 
grades slowly with increasing band indices, since the en- 
velope functions become more and more oscillatory, re- 
flecting the not-plane- wave like nature of Bloch modes. 
It is worth noting, that the Wannier functions optimized 
with respect to the IM and SM criteria are in general not 
equal, even if the global minimum of the SM or the global 
maximum of the IM has been reached. For example for 
the system of dielectric rods in air (TM polarization) the 
SM- and IM-optimized Wannier functions coincide for 
the 1st and the 2nd band (not shown), but not for the 
third band (top right Wannier function in figures [2] and 

BJ. 

Sometimes, it is beneficial to use the time-reversal sym- 
metry of wave operators Ce and Ch to construct real 
Wannier functions. This is possible, since the envelope 
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FIG. 4: Locality (X n ) of the IM-optimized Wannier functions 
as a function of GA generation for the Sq-D structure, TM 
polarization, 3rd band (top) and for the Sq-A structure, TE 
polarization, 5th band (bottom). The dashed green line shows 
the locality of the best-localized Wannier function in each 
generation. Every 100 generations these Wannier functions 
served as a starting point for the subsequent CG optimization 
step (red crosses). On the right-hand side, the modulus square 
of the SM-optimized Wannier functions is shown for an early 
(20000th) and a later (60000) generation, respectively. 



functions of the Bloch modes transform as 
Un-k(r) = <k( r ) 



(12) 



under time reversal (inversion of the reciprocal space), 
and the Bloch phases can be constrained by the condition 



(13) 



By restricting the parameter space of the optimization 
problem in such a way, we found that both SM and IM 
optimization criteria possess multiple extrema, making 
the use of local optimization method not efficient with- 
out special choice of the initial Bloch phases. However, 
using real-valued Wannier functions is of advantage for 
the considerations of the next section. 



IV. CHOICE OF INITIAL CONDITIONS 

Here we propose an analytical expression for a generic 
set of Bloch phases to be used as a starting point for 
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FIG. 5: Modulus square of the maximally localized Wannier 
functions (with respect to the IM) for the Sq-D structure (TM 
polarization). The Wannier center was chosen as 'on-site' for 
the 1st and 5th band and as 'between' for the 2nd and 4th 
band. 
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FIG. 6: Modulus square of the maximally localized Wannier 
functions (with respect to the IM) for Sq-A structure (TE 
polarization). The Wannier center was chosen as 'on-site' for 
the 1st band and as 'between' for the 2nd, 4th and 5th band. 



the optimization procedure. It is based on the following 
theorem. 

Theorem: Suppose that (i) the Wannier functions are 
real- valued, i.e., n k = — n -k, and the Bloch functions 
B n k(r) = e 1( ^ k B n k(r) conform to the following condi- 
tions: Re(B n k(r)) has the same sign (i) for all points r 
in the unit cell and (ii) for all wave vectors k in the first 
Brillouin zone (for each component of the vector Bloch 
function) . Then maximizing the IM of the Wannier func- 
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FIG. 7: Modulus square of Bloch-criterion fEq. U5p optimized 
Wannier functions for a square lattice of dielectric rods in 
air (TM polarization). For n=l,2,4 the Bloch criterion is 
equivalent to the IM-criterion for real Wannier functions. For 
n=3,5,6 the Wannier functions show at least a tendency to 
localize around the Wannier center. 



tions, 



d 2 r|W nR (r)r 



uc 



(14) 



with respect to the Bloch phases n k is equivalent to 
maximizing the IM of the real part of the Bloch functions 
in the first unit cell around the Wannier center, 



X* ({<A„k}) = / d 2 r(Re(e i ^B„ k ))^ 
JVC v 7 



(15) 



for each wave vector in the first Brillouin zone separately. 
This theorem holds in any spatial dimension, even though 
in this paper we consider two-dimensional systems only. 
Note that in the above expressions we have set the weight 
factor X(r) appearing in Eq. ([TO]) equal to unity for sim- 
plicity. 

The proof is straightforward. Due to the translation 
property of the Wannier function, (|8j), it is sufficient to 
prove the equivalence for R = only. In this case 



keBZ 



W n0 (r) 

It is easy to recognize that the equality 

/ d 2 r|Re(W„ (r))| = ^= V f d 2 r|Re(B„ k )| 
Jvc V-N kgBZ Jvc 



(16) 



holds, if the conditions (ii) and (iii) are fulfilled. Since 
the Wannier functions are chosen to be real- valued (i), 
the functional on the left-hand side is maximized by the 
same set of </> n k as Z^({0 n k}), which proves the theorem. 



The maximization problem JTSJ) can be solved analyt- 
ically leading to the following set of Bloch phases 



tan(2(/>k) = 



where B n ^ 



- J uc ^ 2 r2Re(ff nk )Im(ff nk ) 
/ uc d 2 r{Re(5 nk ) 2 -Im(5 nk ) 2 }' 



(17) 



1/2 



B n k(r) • B n k(r) is the (complex) am- 
plitude of the vector Bloch function. Relation ([TT]) de- 
fines the Bloch optimization criterion. Note that the 
phases which are transformed by 0k — » </>k + tt, or equiva- 
lently e i</)k -> e [ ^+^ = -e^ k , also fulfill ([T7J). Further- 
more, the phases which fulfill ([TT]) also obey ([T3]) and, 
therefore, the Wannier functions optimized with respect 
to the Bloch criterion JTSJ) are real-valued, as initially 
assumed. 

Due to the strongly oscillatory nature of Bloch func- 
tions in coordinate space (especially for the higher order 
bands), condition (ii) will hardly be fulfilled exactly in 
a realistic system. But the conditions (i) and (iii) can 
always be fulfilled simply by correct choice of the phase 
factor sign. In figure [71 an example of the Bloch criterion 
optimization is shown for the Sq-D structure (TM polar- 
ization). For the 1st, 2nd and 4th bands Bloch functions 
fulfill Bloch criterion conditions at least approximately, 
and an analytical set of Bloch phases leads to well local- 
ized Wannier functions. Even in the case when the Bloch 
functions do not maximize the Bloch criterion Eq. ([T5]) 
exactly, the Wannier functions optimized with respect 
to ([T5|) exhibit a tendency to localize around its Wannier 
center. The same tendency has been obtained for all four 
test structures in both fundamental polarizations. This 
suggests to use the analytical set of Bloch phases ([TTj) 
as an initial set of phases for the numerical optimiza- 
tion, using either second moment or integrated modulus 
methods. This should help avoiding the local minima 
trapping problem and can reduce the computation time 
considerably. 



V. WANNIER FUNCTION QUALITY 

In this section the quality of IM optimized Wannier 
functions is demonstrated by using them as a basis set 
for photonic crystal defect structure analysis. In fig- 
ure [8] the reconstructed tight-binding band structure of a 
square lattice of dielectric rods in air (TM polarization) 
is shown. The number of lattice sites taken into account 
in the nearest neighbor approximations were increased 
successively. The deviation of the reconstructed band 
structure from the original band structure decreases by 
increasing the number of next neighbors taken into ac- 
count. By restricting to lattice sites separated by up to 
four lattice constants the band structure is reproduced 
well except for small deviation at higher bands and sym- 
metry points. The slight mismatch for higher bands is 
due to the fact that the higher band Wannier functions 
are not as well localized as the lower ones. 
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FIG. 8: Band structure of a square lattice of dielectric rods 
in air (TM polarization). The dashed curve is the original 
band structure, the solid dots are the reproduced band struc- 
ture obtained by within Wannier function formalism. On the 
right side the set of next neighboring sites surrounding Wan- 
nier center R = is shown for different nearest neighbor 
approximations. 



We further calculated the modes and frequencies of a 
point defect structures which consist of a single rod with 
deviating permittivity €d e f at Udef m the Sq-D struc- 
ture. The light propagation inside the band gap between 
the first and the second band is forbidden and the forma- 
tion of a localized defect mode is possible. We split up 
the total permittivity e(r) = e p (r) + 5e(r), into a periodic 
term e p (r) which corresponds to the unperturbed system, 
and a defect term c5e(r) = (cd e f — e p )6(r — |Rd e / — r|). 
An expansion of the z-component of the electric field in 
terms of maximally localized Wannier functions leads to a 
generalized eigenvalue problem with sparse matrices [HI . 
The defect mode frequencies are obtained from the solu- 
tion of this sparse system (Fig. [9]). The first eight Wan- 
nier functions have been taken into account. There are 
monopole like defect modes for a defect rod permittivity 
of €def < 12 and doubly degenerated, dipole like modes 
for higher permittivities in the defect rod €d e f > 12. The 
results are in complete agreement with plane wave cal- 
culations (23) . To analyze contribution of the individual 
Wannier functions to the defect mode, the band index 
contribution 



R 



VI. CONCLUSION 

The procedure to construct maximally localized Wan- 
nier functions by Bloch phase optimization was analyzed 
for several two dimensional photonic crystals for both 
fundamental polarizations, using two different locality 
measures. Although the stochastic, genetic algorithm is 
numerically too costly for routine application, it has, as 
a global optimization method, provided us an important 
benchmark to judge under which conditions the faster 
and less memory intensive, but local conjugate gradient 
method finds the global optimum of a given locality mea- 
sure. We found that the commonly used second moment 
locality measure has generically multiple extrema, which 
makes it difficult to construct maximally localized Wan- 
nier functions by local optimization techniques. One may 
conjecture that this multiplicity results from the complex 
oscillatory behavior of the (not yet optimized) Wannier 
functions at large distances from the Wannier center, 
which makes the dominant contribution to the second 
moment. This led us to propose a new locality mea- 
sure which is controlled by the behavior close the Wan- 
nier center, the integrated modulus square measure. We 
showed numerically by comparison of conjugate gradi- 
ent and genetic algorithm optimization that this measure 
does not feature multiple extrema and is, therefore, SUit- 
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where the normalization factor is given by M = 
^2nn\-^nn\ 2 is shown in figure [TUI for the defect modes 
€def = 1 and Cdef = 30. One can see, that C n rapidly 
decrease for higher band indices. Since the defect mode 
frequencies are in between the 1st and the 2nd band, only 
the lower band Wannier functions contribute to the de- 
fect modes. Therefore it is justified to cut of the band 
index which reduces the numerical load. 



FIG. 9: Frequencies of the modes in a point defect consisting 
of a single rod with differing permittivity €d e f in a a square 
lattice of dielectric rods in air (TM polarization). The dots 
indicate the results of the Wannier function approach by tak- 
ing the first eight bands into account. They are in complete 
agreement with plane wave calculations (red line) [23|] . At the 
bottom the real part of two defect modes with 6d e f — 1 and 
Cdef = 30 is shown. 
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n 

FIG. 10: The band index contribution C n for the defect modes 
€def = 1 and €def = 30. Since the defect frequency lies inside 
the band gap between the 1st and the 2nd band, only the 
lower band Wannier functions contribute to the defect mode. 

able for fast and efficient local optimization techniques, 
like the standard conjugate gradient method. Because 
this result presumably originates from the local nature 



of the integrated modulus measure, it should hold gen- 
erally, not only for two-dimensional systems, but also in 
three dimensions for photonic as well as electronic lat- 
tices. 

We also presented and tested an analytical formula for 
the set of Bloch phases to be used as a starting point of 
the optimization process. This initial set of Bloch phases 
is suggested because, albeit it does not solve the opti- 
mization problem in general, it does generate maximally 
localized Wannier functions in special cases where the op- 
timization problem can be solved analytically. We expect 
that these two main results may significantly increase the 
efficiency of the Wannier function approach for the de- 
scription of defect structures in photonic lattices. 

This work was supported in part by the Deutsche 
Forschungsgemeinschaft (FOR 557). 
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